Variable-range-hopping in two-dimensional system. 
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Computer modeling of the VRH conductivity in the two-dimensional system has been done 
by kinetic Monte Carlo method, which includes some new elements. Study of the temperature 
dependence of the conductivity, testing of the different scaling relations, and study of the size effect 
show that the Efros- Shklovskii mechanism of the VRH is valid in the slightest details. It has been 
also shown that simultaneous transitions of many electrons are not important. The reasons of 
disagreement with previous works are thoroughly analyzed. 



The concept of the variable range hopping (VRH) of 
localized electrons belongs to Sir Mott jjj. He has con- 
sidered the phonon assisted tunneling of localized elec- 
trons between different sites. Mott has found that the 
typical hopping length for a two-dimensional (2D) case 
R M w (ci/3)(Tm/T) 1 / 3 increases with the decreasing 
temperature. This is the origin of the term "VRH" . Here 
Tm = 0M(go a2 )~ 1 , a is the localization length of an elec- 
tron, go is the density of states at the Fermi level, (3m is a 
numerical coefficient. The conductivity obeys the Mott's 
law §0 



o~m 



( 7 e 2 /T)exp (-(Tm/T) 1 ' 3 



(1) 



where 7 is a typical phonon frequency of the order of 
KPsecT 1 . 

Electron-electron interaction has added a new chapter 
to the theory of VRH. The study of the interaction in lo- 
calized regime has been initiated by PollakH] and Srini- 
vasan jE). Later on Efros and Shklovskii [[jhave shown 
that the single particle density of states (DS) in 2D case 
linearly tends to zero as energy tends to the Fermi energy. 
This phenomenon, which is called the Coulomb gap(CG), 
is due to the long range part of the Coulomb interaction 
which, in a sense, remains non-screened. At large disor- 
der the DS G(e) = 2\e\/ire 4: has a universal form. The 
DS is given by the only combination of the energy and 
electron charge which has a proper dimensionality. In 
fact, the CG results from the Coulomb law and from the 
discrete nature of electron charge. 

Simple arguments based upon single-electron excita- 
tions show that the VRH with the Coulomb interaction 
(CI) obeys Efros-Shklovskii (ES) law @ 



( 7e 2 /T)exp (-(To/T) 1 / 2 



(2) 



where To = Poe 2 /a. The percolation approach (See ref- 
erences in ||) gives flu = 13.8, ^0 = 6.5. The hopping 
length R c m (a/4) (To/T) 1 /*. 

Many experimental works on the 2D structures report 
the ES law, but very often the data are ambiguous be- 
cause one should have a really large interval of a to dis- 
tinguish between the two laws. We think that the most 
important results have been obtained recently by groups 
of Jiang and Dahm and Adkins B. Using an artificial 



screening provided by metallic electrode, parallel to a 
plane with the 2D electrons, they have proved that the 
VRH transport reflects the crucial feature of the CG, the 
sensitivity to the long range interaction. 

The theoretical understanding of the VRH is contro- 
versial. Efros and Shklovskii argue that in 2D case only 
the single electron transition with the typical hopping 
length are important ||[)| and the VRH has the same uni- 
versal nature as the CG. Pollak and his followers (h] [r| 
claim that both sequential and simultaneous transitions 
of different groups of electrons are very important for 
the VRH, so that the physics of the CG which is valid 
for the long hops only is not relevant. These claims are 
based upon the theoretical idea that simultaneous tran- 
sition of many electrons may provide a gain in the acti- 
vation Jul . Computer simulations of the VRH and the 
similar phenomena made by this group also show impor- 
tance of the many-electron transport both in 3D and 2D 
cases but no new T-dependence of the conductivity has 
been proposed. The simulation of our group |3j reveals 
scaling typical for the CG, but in this work a simplified 
model has been used which does not take into account 
any simultaneous transitions and underestimates sequen- 
tial transitions with the short length. 

We present here the results of computations which 
takes into account all the sequential transitions and the 
most important simultaneous transitions of two electrons. 
We show that the VRH both in the CI case and without 
interaction (M-case), closely obeys the laws Eqs. 
and different scaling laws which reflect the physics of the 
both cases. We also confirm earlier prediction that simul- 
taneous transitions of many electrons are not important 
for the VRH in 2D case at low temperatures and analyze 
the discrepancy with the simulations by Pollak group. 

We consider the standard lattice model with the 
Hamiltonian 



H = ^ <pin t + - ^2 — ( n i ~ v ) i n j - v) j 



(3) 



where = 0, 1 are the occupation numbers. The 
quenched random site energies fa are distributed uni- 
formly within the interval [—1,1], v is the average occu- 
pation number, which is taken to be 1/2 everywhere be- 
low. The magnitude of the quenched disorder is enough 
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to provide the universal CG at all energies important in 
our temperature range Here and below the lattice 
constant is a length unit and the nearest Coulomb en- 
ergy which is equal to the amplitude of the disorder is 
both energy and temperature unit. In the M-case the 
interaction term in Eq. (||) is neglected. 




FIG. 1. The temperature dependence of the conductivity 
as obtained by simulations for the M-case (a) and for the 
Cl-casc (■). Averaging over different disorder realizations is 
performed to get the error bars equivalent to the sizes of the 
symbols. Dashed line shows the fit by the Mott's law Eq. (|l|), 
while solid line represents ES-law Eq. (Q). The localization 
length o = l. The inset shows the test of the scaling laws Eq. 
(§). Here Ta 2 = 0.04 for the M-case (a) and Ta = 0.1 for the 
Cl-case (■) (See the explanations in the text) 

To simulate the VRH by the kinetic Monte Carlo (MC) 
method we add an electric field E and consider the pe- 
riodic boundary conditions (See details in Ref ||). To 
get conductivity we calculate the dipole moment due 
to electron transitions in the -E-direction and divide it 
by the number of MC steps and by E. This assumes 
that our time unit is I/7. We check that the result is 
_E-independent. In fact we take E from the condition 
ER c ,m = 0.1T. The array size is taken 100 x 100 for the 
temperatures above 0.05 and 200 x 200 for the temper- 
atures below 0.05. In the first approximation we neglect 
simultaneous transitions of many electrons. The new el- 
ement of our algorithm is that we do not include tunnel- 
ing probability into MC process |lq] . Instead we choose a 
pair of sites i and j with the probability exp(— 2r*y /a) jZ\ , 
where Z\ is the normalization factor. Then the transition 
is rejected if both sites are occupied or empty, and finally 
it is performed with the probability 1/(1 + exp(— e^/T)), 
where is the energy difference between the two con- 
figuration. This new scheme dramatically decreases the 
runtime of the simulation and leads to L 2 , rather than 
L 4 , algorithm in number of MC steps. Note that for the 



Cl-case each MC step consists of L 2 operation itself due 
to recalculations of site energies after each transition. 

To reveal and compare the physics of the VRH in both 
CI- and M-case we perform both simulations. Fig. 1 
shows that the laws Eqs. @,(0) are reproduced very 
good for both CI- and M-cases. From solid and dashed 
line one gets (3q = 5.8 and (3m = 13.0 respectively, that 
agrees well with the percolation values cited above. Note 
that the temperature where two curves become close to 
each other marks the upper temperature limit of the 
VRH. At even higher temperatures the two curves co- 
incide and the conductivity decreases with the increas- 
ing temperature which manifests transition to a diffusion 
regime. 

To test the origin of the VRH we check the scaling 
laws which reflect an intimate physics of both CI- and 
M-cases: T and Tm respectively are the only character- 
istic temperatures (energies) in the problem. These laws 
can be written in the universal form 

<?m = (je 2 /T M )fM (T M /T) , , , 

ctci = he 2 /T )fci (T /T) , W 

where Jm{x) and fciix) are some functions. It follows 
from the scaling laws Eq. ^ and the above definitions 
of Tm and Tq that <tji//o 2 and &ci/ a are functions of 
Ta 2 and Ta, respectively. Thus, if one changes both 
T and a in the M-case in such a way that Ta 2 is con- 
stant then crjv//a 2 will be also constant. For Cl-case one 
should change T and a keeping Ta constant to get aci / a 
constant. This test has been successfully performed and 
the convincing results are shown in the inset of Fig.l. 
Note that the small deviation of a /a 2 from the constant 
value for the M-case is due to the lattice effect. Since 
Rm < Rc, the lattice effect in the M-case appears at 
larger a than in the Cl-case. 

In Fig. 2 we present the results for the size effect in the 
Cl-case. We have calculated the average value of the con- 
ductivity av= (J/, /(Too, where 07, is the logarithmic aver- 
age of the conductivity over arrays with different disorder 
at given L and T, cr^ is the conductivity at the maxi- 
mum L at a given T. We have also calculated dispersion 
dis= Sa/a, where 6a is the dispersion of the conductiv- 
ity and a is the arithmetic average of the conductivity at 
given L and T. 

We did these computations not only to get an idea how 
large the size of an array should be to obtain the macro- 
scopic value of the average conductivity and what value 
of the dispersion one should expect at a given L and T. 
Another reason is that the L- and T-dependences of the 
average and the dispersion give us a new insight into the 
physics of the VRH. One can see from Fig. 2 that the 
dispersion is a function of only one argument TL, as it 
follows from the ES theory. Indeed, the correlation length 
L c of the percolation network is |tf| L c = Rc{To/Tf/ 2 . 
Since the exponent of the correlation length v = 1.33, 
one gets that Lc is very close to 1/T. Thus, L/Lc ~ LT 
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which explains the L- and T-dependence of the disper- 
sion. 

The behavior of the average is completely different, ft 
does not change until very low values of TL. At even 
smaller TL it decreases dramatically and it is not any 
more a function of TL. We believe that it is a func- 
tion of L/Rci ~ Ly/T, and the decrease of the average 
conductivity appears due to the hard gap in the DS. It 
is known Jl7[| that the hard gap between occupied and 
empty states results from the size effect and its value is 
of the order of 1/L. At L < Rc the energy stripe of the 
VRH is completely within this gap. We suggest that this 
is the reason of the decrease of the average conductivity. 
Fig. 2 confirms qualitatively this suggestion. Unfortu- 
nately, we are unable to check it quantitatively. 




FIG. 2. The size dependence of the dispersion (left axis, 
opened symbols) and of the logarithmic average of the con- 
ductivity (right axis, solid symbols) at T — 0.1 (A, A) and 
T — 0.05 (□,■) for the Cl-case. Inset: the conductivity at 
the presence of the simultaneous transitions of two electrons 
for the Cl-case. The contribution to the conductivity due to 
single hops (•) and due to double hops (o) at a = 0.5 should 
be compared. The conductivity with single hops only at a — 1 
(□) is also shown. 

Now we come to the simultaneous transitions of many 
electrons in the interacting system. We consider only 
two-electron transitions, which are the most important. 

Suppose that in the initial state electrons occupy sites 
1 and 3 and in the final state they are in the sites 2 and 4. 
The activation energy of this process £1234 is the energy 
difference between the final and the initial states. Due 
to Pollak it may be smaller than the sum of the energies 
of the single-electron transitions. Moreover, the increase 
in the activation exponent may beat an extra tunneling 
exponent, which appears because of the tunneling of two 
electrons. That is why Pollak and his followers insist on 
many-electron processes. 

First we calculate the probability of the process. Let us 
consider the case when ryi , <C R, where R is the dis- 



tance between pairs (1,2) and (3,4). Following Ref. (h]] 
we concentrate on the single phonon transitions. It is 
obvious that these transitions only cannot provide any- 
thing but independent hops of electrons. The result is 
different if one considers the CI as well. We think that 
the amplitude of the most important process is described 
by the second order perturbation theory. It consists of 
a product of two matrix elements. The first one gives 
the electron transition in one of the close pairs, say (1,2), 
due to the phonon absorption without energy conserva- 
tion. This element has a factor exp(— r^/a). The system 
comes to the intermediate state, where sites 2 and 3 are 
occupied. The second element describes the transition of 
the second electron from 3 to 4 in the external electro- 
static potential created by the dipole (1,2). This element 
has a factor exp(— r^/a)/ R 3 . The independence is due 
to the dipole-dipole interaction. Use of the Golden rule 
gives the correct energy conservation law and a factor 
~ exp(— 2(ri2 + r34)/a)/R 6 . The result is the same as 
for the Forster diffusion of molecular excitons Jf^] , very 
similar but not identical problem. The transition rate 
used in Ref. [|fl]-[l3| and many others is taken to be R- 
independent at large R. We cannot understand this re- 
sult because to get the common energy conservation law 
one should assume that the remote pairs can "talk" to 
each other, which means there should be an interaction 
between them with some independence. We argue be- 
low that the contradiction between our results originates 
from the independence. 

To accommodate this fast decay with R we assume that 
the simultaneous transitions occurs only if the initial or 
final position of the second electron is within the spheres 
with the radii Rc from the initial or final position of the 
first electron. Note that this assumption overestimates 
the role of the double electron transitions in the most 
important case when one electron makes a long hop Rc 
and the other one makes a short hop. The interaction be- 
tween such pairs decreases (as 1/i? 4 ) at distances smaller 
than Rc- 

To simulate the double hop transport we use the same 
idea as for single hop case. This time, however, we 
choose a pair of sites i and j for single hop transi- 
tion with the probability p± — exp(— 2rij / a) / Z-i, and 
sites i, j, k and I for double hop with the probability 
P2 = exp(— 2(ry + t^i) / a) / Z2, where Zi is found from 
the normalization condition for the sum of two probabil- 
ities. 

If the double hop is chosen, we check the occupation 
for both pairs of sites and (k,l) and then use the 

same Boltzman factor as before with eijki instead of e^. 
The algorithm is also L 2 in the number of MC steps. 
However, since the number of possible transitions of two 
electrons is much larger than the number of single tran- 
sitions, one can show that there are w 4i?^ exp(— 2/a) 
attempts of double transition per one attempt of a single 
one, which makes simulation much more time consuming. 
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The contribution of single and double hops to the total 
conductivity can be separated by collecting their dipole 
moments in the different cells. We have checked that the 
single hop contribution does not change in the presence 
of double hops. The results of the simulation are shown 
on the inset of Fig. 2. We are able to calculate double 
hop contribution for the case a = 0.5, because at larger 
a the value of Rc is too large. The single hop contribu- 
tions are shown for both a = 1,0.5. Note that a small 
deviation from the universal behavior at high tempera- 
tures is a result of the lattice effect for a = 0.5. One can 
see that the double hop contribution is from one to two 
orders of magnitude smaller than the single hop contri- 
bution at VRH region of temperatures. Moreover, the 
relative fraction of double hops in the total conductiv- 
ity decreases with decreasing temperature. This result 
does not support Pollak's prediction [jl0| of the gain in 
the activation energy for simultaneous transition of two 
electrons as compared with the single electron transition. 

The computational data, obtained in the papers 
0-H and many other works of this group show the 
importance of the many-electron hops which contradicts 
to our result. In our opinion the i?-independent proba- 
bility of the simultaneous hops assumed by this group is 
responsible for the contradiction. 

Indeed, if in a large array each pair has an equal prob- 
ability to make a transition with any other L 2 soft pairs 
in the system the number of transitions for the pair per 
time is proportional to L 2 . This is a giant overestimate of 
the role of many-electron transitions. The abnormal size 
effect should be a manifestation of this mistake. Unfortu- 
nately, we have found only one work p2| , where authors 
take care about the size effect and this work deals with 
the 3D case. However, we think that the infinite range 
interaction between pairs should provide the same effect 
in the 3D-case. Looking at Fig. 3 of Ref. one can see 
an exponential size effect without any signs of a satura- 
tion. The resistivity increases with the size because the 
algorithm include the selection of states with the highest 
transition rate among them. The selected states must 
consist of soft pairs which are far from each other and do 
not contribute to the dc conductivity. 

Another reason for contradictions is that the method 
used in all these papers is applicable to the mesoscopic 
systems only. It consists of selecting a number of low- 
energy states of the total system to study transitions be- 
tween them. It looks nice because it takes into account 
all many electron transitions, while we could afford the 
simultaneous transitions of two electrons only. The prob- 
lem is in the necessary number of states, which provides 
reasonable thermodynamic and kinetic description. This 
number increases exponentially with the number of sites 
in the array. Due to the sensitivity of the CG to the long 
range interaction the size of the array L, which provides 



macroscopic regime, increases with decreasing T. One 
can estimate from Fig. 2 that L w 3/T gives a reason- 
able dispersion. For example, T= 1/400 is the middle of 
the temperature interval in Ref. jflj. It follows that L 
should be at least 1000. Instead, authors compute few 
samples with different disorder at L » 22. These samples 
have very different values of conductivities and different 
patterns of the T-dependences (See Fig.l of Ref. Jil"|). 
This is a typical mesoscopic result. 

Finally, our computer simulation confirms the ES the- 
ory of the VRH transport in 2D systems and shows that 
the attempts of its revision have no ground. 

We are grateful to O. Biham and B. I. Shklovskii for 
the interest to this work and to E. I. Rashba for thorough 
discussion of the theory of simultaneous transitions. The 
work have been funded by the US-Israel Binational Sci- 
ence Foundation Grant 9800097, the computations have 
been made in CHPC of the University of Utah. 
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